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The concept of stochastic matrix product states is introduced and a natural form for the states 
is derived. This allows to define the analogue of Schmidt coefficients for steady states of non- 
equilibrium stochastic processes. We discuss a new measure for correlations which is analogous 
to the entanglement entropy, the entropy cost Sc, and show that this measure quantifies the bond 
dimension needed to represent a steady state as a matrix product state. We illustrate these concepts 
on the hand of the asymmetric exclusion process. 
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This paper explores parallels between the many-body 
description of non-equilibrium steady states of classi- 
cal stochastic processes and ground states of strongly 
correlated quantum many-body systems. Classical non- 
equilibrium steady states are typically much richer than 
their equilibrium counterparts and can exhibit interest- 
ing behavior such as the presence of a current, non- 
equilibrium phase transitions and entire phases with a 
diverging correlation length pH f2l. There are strong 
parallels between the many-body description of non- 
equilibrium classical stochastic spin systems [3 , and it 
turns out that matrix product states (MPS) play a very 
important role in both fields [11 15]. 

The theory of entanglement, developed in the context 
of quantum information theory, has recently proven to 
yield valuable new insights in the nature of the wavefunc- 
tions arising in such strongly correlated quantum many- 
body systems. One of its main virtues is that it allows to 
quantify the amount of quantum correlations in a many- 
body system in terms of so-called area laws [6]-[8] and re- 
late this to the structure of the underlying wavefunctions 
representing the ground states of quantum spin Hamilto- 
nians [9l[T0]. The ground states of quantum spin chains is 
well described by the class of MPS, since for local Hamil- 
tonians the associated Schmidt coefficients of the ground 
state decay very fast. A MPS approximation with finite 
bond dimension can be obtained by setting the smallest 
Schmidt coefficients aj equal to zero. This approxima- 
tion is justified whenever an area law is satisfied [9 . If 
we define pi = af, then those pi form a probability dis- 
tribution. The entropy of this distribution quantifies the 
dimension of the matrices that has to be chosen in the 
MPS approximation such as to obtain a given fidelity. 
For such pure states, this entropy has multiple opera- 
tional meanings, ranging from the amount of distillable 
Bell states by local operations to the maximum amount 
of secret information that can be sent in a cryptographic 
setting [TT] . For mixed states and in the presence of clas- 
sical correlations the situation is more complicated. A lot 
of research in quantum information has concentrated on 
resolving the relationships between classical and quan- 
tum correlations in mixed quantum states. A measure 



of particular importance to our approach is the entan- 
glement of purification It quantifies both quantum 
and classical correlations in terms of the number of maxi- 
mally entangled states needed to asymptotically generate 
the state. 

The main topic of this paper is to explore the classi- 
cal non-equilibrium analogous of the quantum notions of 
entanglement entropy, area laws and the density matrix 
renormalization group which justify the use of MPS in 
the quantum setting. The main technical difficulty arises 
from the fact that classical probability distributions are 
normalized in the Li norm \pi\ = 1), while quantum 
states are normalized in the L2 norm \^|^^\'^ = 1). 
This can partly be overcome by working with a subclass 
of MPS where all matrices only contain non-negative 
entries; we call this subclass stochastic matrix product 
states (sMPS). The concept of mutual information, de- 
fined for classical bipartite distributions pab where A 
and B will represent the variables or spins on both halves 
of a chain, plays a role analogous to the entanglement 

entropy: I{A : B) = Eab Pab logs (^^)- 
diately gives an upper bound to the error made when 
approximating pab by a product of its marginals, since 
Wpab - PaPbWi < 21n(2)/(A : B) [13. Just as in the 
quantum case, one would expect that the global non- 
equilibrium steady state probability distribution of the 
stochastic process can be represented as a stochastic ma- 
trix product state (sMPS) with small bond dimension if 
this mutual information is small. However, more subtle 
measures are needed in the case of stochastic processes, 
and we will introduce the notion of entropy cost to quan- 
tify the bond dimension needed for the corresponding 
sMPS. 

Definition 1 The entropy cost Sc for a bipartite prob- 
ability distribution P(x, y) is given by: 

Sc= min S({px}) (1) 

where S{{px}) = — X1;^Pa log2 (Pa) ^-^ the Shannon 
information of {px}, and where the optimization is 
over all probability distributions px and over all con- 
ditional probabilities Pa and Pb for which P(x, y) = 



2 



j:^PA{x\X)PB{y\\)Px. 

The entropy cost bears a lot of resemblances to the no- 
tion of common information introduced by Wyner in the 
context of cryptography and classical information the- 
ory [M], and this entropy cost serves as an upper bound 
to the common information. The entropy cost can be 
thought of as the classical analogue of the entanglement 
of purification, and the probability distribution px plays 
a role analogous to the Schmidt coefficients in the quan- 
tum case. 

To give an example, consider the equilibrium system 
defined by the classical Ising model H = YliLi^ 
with Si = ±1. The equilibrium distribution p{{si})i = 
ex.p{—l3H)/Z can be written in terms of a very simple 
MPS with D = 2 [15]; the entropy cost Sc can then 
easily be calculated exactly in the case N = 2 and is 
given by 



Sc = - (e"^ cosh(/3) log2 (e-^ cosh(/3)) + 
e-^ sinh(/3) log2 (e-^ sinh(/3))) 



(2) 



As expected, this function monotonously increases from 
to 1, i.e. from the paramagnetic without correlations 
to the ferromagnetic region with strong correlations. A 
different measure that was recently investigated is the 
shared information [16 . It has been shown that it obeys 
an area law for several non-critical stochastic models and 
that critical behavior can be identified by logarithmic 
corrections. 

Let us next define a D-dimensional sMPS describing 
a classical probability distribution of A/'-spins each of di- 
mension d] obviously, those sMPS were already exten- 
sively used in the literature, and we are just formalizing 
the definition here. 

Definition 2 A stochastic matrix product state 
(sMPS) is given by: 



\Pd) 



^ {L\Bl...Bl\R)\h...i,,), (3) 



where we only consider real matrices that are x Z^^+i 
dimensional, with Dk < D, and additionally fulfill the 
requirement [^f^]^ > for every element individually. 
This ensures that all the weights of the distribution are 
positive after contraction. The left and right vector (L | 
and I R) are also elementwise positive and can be ab- 
sorbed into the matrices B}^ and Bf^^ , which corresponds 
to choosing Di = Dn = ^- Furthermore, we require \pd) 
to be normalized in the Li norm, || \pd) ||i = 1- 
Every multipartite probability distribution of a chain of 
discrete variables can obviously be written in the form 
(|3| if we allow for a sufficiently large matrix dimension 
2)max^ i.e. exponential in the number of sites. 

Let us next show how the entropy cost can be 
calculated. Upon inserting a partition of unity 
I = Ea=iI^)('^I in ([I), we can write \p) = 



EAEK}(^l^^..^fJA) (A|5f+i^...<|i?)|{z4). 

Now k N 

PA=(iinc'"'iA)(Ai n (4) 

n=l n=/c+l 

defines a new probability distribution if we chose C^^^ = 
^i^ B^l^^ (i.e. the transfer matrix). This allows us to 
rewrite the MPS as 

\P) = Y. PA{{in}neA\X)PxPB{{in}neB\\) (5) 



A=l ii 



where 



PA{{in}neA\>) = 



PB{{in}neB\>^) 



{L\Bl...Bl\\) 

{L\Y{Uc\n\x) 



(6) 



The probability distribution {p\\ sums up to one due to 
the normalization of (|3|, and PA({^n}|A) and PB({^n}|A) 
can be interpreted as information channels. Note that 
there are several partitions of unity that will give rise to 
valid Pa^Pb and px^ and therefore the decomposition ([5| 
is not unique. It is probably a NP-hard problem to find 
the optimal decomposition that minimizes the entropy 
cost, and in practice we will therefore rely on the con- 
struction that was just given for finding upper bounds 
to it. It is also easy to find lower bounds to the entropy 
cost: 

Lemma 1 For a given distribution \p) the mutual in- 
formation I {A : B) is bounded by the entropy cost Sc- 
I{A:B)<Sc- 

Proof: By virtue of ([5| we can focus on calculating 
the mutual information /(A : /i) of the distribution 
P(A,/i) = p\5\^^^ since the full distribution can be read 
asP^B = EA,M ^A({^A}|A)PB({iB}|/i)P(A,/i). This cor- 
responds to a source that generates two outputs, which 
are then transformed by the channels Pa and Pb- The 
data-processing inequality guarantees that the mu- 
tual information of the processed source I{A : B) = 
I{Pa{P) ' Pb{P)) ^ ^(A : /i) is smaller than the mu- 
tual information of the source itself which is equal to 
its entropy /(A : ji) = S{{px}). The decomposition 
[Pa({z^}|A), Pe({i^}|A),p(A)] is not unique, therefore 
the bound is improved by taking the minimum over all 
decompositions. □ 

The decomposition of the sMPS as given in ([5| suggests 
the existence of the following (non-unique) normal form: 

\p)= J: A^pmA^...P^--^A^\n...^r,). (7) 

ii...iiv 

Here the matrices P^^^ represent diagonal matrices with 
probabilities {pa^} sorted in decreasing order. The ma- 
trices = T.i^A^^\ where [A^^^]] > 0, combined 
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with the p[^]'s behave as stochastic matrices. That is, 
pk-i^lk] ^ Qkpk ^ gT^ ^Yiem and Sk de- 

note different stochastic matrices. 

To see that every MPS distribution can be written this 
way consider the fohowing scenario: We start by intro- 
ducing the first bipartitioning between the first two sites. 
After the necessary normalization we proceed to the next 
site and perform the same procedure renormalizing the 
resulting matrices by the total contraction of the two 
halves of the chain. Proceeding along the chain results 
in the desired form. 








ik+i 




ik+2 ^ 



(g) 



1Z 



FIG. 1: Pictorial representation of the natural sMPS decom- 
position: Image (a) can be seen as the graphical representa- 
tion of eqn. From a given source Pa, the correlations are 
distributed via the two channels on the left and right. The 
normalizing factor is included in the A's. In (b) the analogy 
to the quantum MPS becomes evident for the decomposition 
as given in eqn. ([7]).The probabilities in the matrices P^^^ 
are the analogous of the singular values which arise upon a 
Schmidt decomposition of the quantum state {17\ . 

This representation ^ enables us to give a good estimate 
on the error measured in the Li norm which is made 
upon truncating the dimension of the source space, i.e. 
neglecting probabilities smaller than a given value along 
each bipartition: 

Lemma 2 For every multipartite distribution \p) there 
exists a MPS \pd) of the form ^ with dimension D, 
such that 
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\Pd) 



N-1 

E 

fe=i 



,i<2 } ^ek{D), 



j-^max tt,-! 
where ek{D) = Ea=d+i^a • 

Proof: We can always write \p) as a distribution of the 
form Q with a D^^^ = d^ . We now introduce another 
MPS pd) natural form with a bond dimension of D. 
Let \pd) = \Pd) /II \Pd) 111' where \p})) is the unnormal- 
ized probability distribution which arises from neglecting 
along each cut all the probabilities {pAfcl^+i • write 
\Ph) = Eii,}Af^P*^'^---AfJ\{h}Li)- Note, that 
if II |P) - \Ph)h < e, then \\\p) - \pn)h < 2e 
due to the triangle inequaUty. Since 1^1^) arose 
by only neglecting positive numbers, we may write 
IIIp) - = cWpi...Cl^l - CWP*1...C[^1 = 

I _ (^iN-l |*)^[N] ^ ||(^;N-l I _ |*)||^^ 

Here we have used the fact that all summands are 



positive and we defined (/^ | = C^-^^P^-^L . . p[^] as 
weh as (/^l* = CWP*W ...P*[^]. The difference 

- <^M*lli = ii(i...i|(p[^]-p*[^Olli = 

Y.al=D+iP>^ is simply given by ei(P>). Note that 
due to ^ = Proceeding to 

calculate the difference for other k we find that 

\\(l[k] I _ (l*[k] 1 < II (^(l[k-l] I _ (i*[k-l] 1) c[k]p[k] + 

||(r[fe-i]|cw(pw-p*w)||^ < Y!rJi^n{D) + 

||(1...1| (Pl*^! -P*!*^])!!^. The last summand cor- 

responds exactly to Y.ak=D+iP(y-k = ^k{D), which 
completes the proof. □ 

We have therefore proven that there exists an efficient 
parameterization of the steady state in terms of a sMPS 
with low bond dimension if there exists a parameteri- 
zation of this steady state for which the entropy cost 
with respect to all bipartite cuts is small: if this is the 
case, then "^^^i ^k{D) can be made small by following 
the arguments outlined in [9j. This is the analogue of 
the quantum case for which the existence of an area law 
implies the existence of an efficient representation in 
terms of MPS. Note however that the classical statement 
is a bit weaker, as the same normal form has to be 
used with respect to all bipartite cuts, and there is no 
guarantee that the same parameterization is optimal for 
all of bipartitions. 

To make the investigations concrete, we consider the 
non-equilibrium steady state of the asymmetric exclusion 
process (ASEP) [18]. This process is modelled by a 
chain of sites labeled k = 1 . . .N occupiable by hardcore 
particles, i.e. classical spins ik G {0,1}. The particles 
are only allowed to hop to the right, and this only if the 
next site is empty. To drive the system, particles at the 
left are injected with a given rate a and removed on the 
right with a rate {3. This process exhibits three phases 
determined by the inflow a and outflow {3. As was shown 
in [18j the steady state \p) of the corresponding master 
equation can be found exactly in terms of a MPS, albeit 
one for which the matrices are infinite dimensional and 
not necessarily positive. 

The unnormalized steady state solution of the ASEP 
master equation ^] is given in MPS form ([3| with 
site independent matrices Pq — • • • — — ^ 
B\ = ... = = V and boundary vectors (P | , | P) , 
where all those satisfy the algebraic constraints 
V£ = £^V, {L\S = ^{L I and V\R) = ^ \ R). Except 
for the special case when a + /3 = 1, the representations 
obeying those constraints are infinite dimensional. As 
the total occupation number of particles is limited by 
we will construct truncated representations of a given 
dimension P = + 1 that still reproduce the exact 
solution for a chain of length N] the entropy cost can 
then immediately be upper bounded by the logarithm 
of the system size, just as in the case of critical quan- 
tum spin chains. These representations are of the form 
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^ = EL=o [Bk, I i) {j \+En=2\n){n\ + \n){n-l\ and 
= EL=o I i) {j \+En=2\n){n\ + \n-l) (n |. 
With the left and right vector {L \ = Yll^=o I and 

I ^) — X]n=o '^n I ^) • Here A and 5 are 2-dimensional 
matrices fulfihing AB + cr~(j+ = A -\- B (here (j+ and 
cr~ are the Pauh raising and lowering operators), and 
the 2-component vectors v and w must be chosen to be 
eigenvectors: Av = {1/P)v^ wB = {l/a)w. 
The construction given in ([5| gives us a way of bounding 
the entropy cost as a function of the parameters a, /3 and 
N. We varied over all A, 5, w satisfying the algebraic 
relations, and found that the minimum is obtained 
for 3 different solutions depending on the parameter 
range of a and f3. The solution is depicted in Fig|2] If 
a + /3 < 1 and /3 < a, then A = | 0) (0 | + 1 1) (1 |, 
B = 1/(1 - /3)|0)(0| + 6|1)(0| + l/a|l)(l| , 
w = {a{l-P) 6)andv^ = (l 0) with b = ^/T^^o^^ . 
However, if a + /3 < 1 but /3 > a, the optimal solu- 
tion can be obtained from the previous one by the 
replacements, A = B^{a ^ f3), B 
"^(q/ ^ For a 4 



(a) 



W = V" 



V = w 



solution is given by A = 1//3 | 0) (0 | 
B = l/a\ 0) (Q| +a|l)(Q| + | 1) (1 
and a = ^^/l/a -\- 1/ (3 — l/aj3. 



- A^ [a ^ and 
/3 > 1, the optimal 
Fa|0) (1| + |1)(1|, 
, = w = (I 0) 





FIG. 2: (color online) Figure (a) The mutual information of a 
chain length N — 20 for different inflow parameters a and (3 
has been calculated numerically. Figure (b) entropy cost for 
a representation with D^^^ — 21. Note the different scales of 
the two plots. The resulting plots clearly reflect the under- 
lying phase diagram of the ASEP [18], both for the mutual 
information and the approximate entropy cost. The entropy 
cost as well as the actual mutual information drop to zero 
along the mean- field line a -\- (3 — 1. The mutual information 
is consistently low throughout the diagram, explaining why 
mean-field approaches have given such good results 1 . 

Conclusions: We have revisited the notion of stochastic 
matrix product states, and showed that a low bond di- 
mension suffices to efficiently parameterize steady states 
of non-equilibrium distribution if the entropy cost in the 
system is low. This opens up the interesting question 
to characterize the conditions under which such steady 
states have a low entropy cost. It would be interest- 
ing to see to what extent this relates to the gap of the 
corresponding stochastic process. This also opens up 
novel ways for constructing numerical renormalization 
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FIG. 3: (color online) Plot (a): Monte Carlo simulations of 
the mutual information for different a and /3 values. The sys- 
tem size was varied from A/" = 2 ... 180. The chain was cut in 
the middle, the length L of the individual blocks are given by 
half of the total system size. The plot shows the exponential of 
the mutual information exp I [A : B). The simulations suggest 
that the mutual information grows only sub-logarithmically 
along the coexistence line a = a + /3<l; Plot (b): The 
logarithm of the Probability distribution {p\} plotted for dif- 
ferent values of a and /3. The dimension of the matrices 8 and 
V is D™'''' = 21, corresponding to a chain of length AT = 20. 
The distributions decay super-exponentially. For all A > 11, 
PA = 0. 

group methods for simulating non-equilibrium systems in 
the line of the MPS algorithms for quantum spin chains 
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